# R Studio version 2024.4.1.748 ("Chocolate Cosmos") on macOS Sequoia (Version 15.5)

library(tidyverse)
library(readxl)
library(dplyr)
library(magrittr)
library(cjpowR)
library(plotrix)
library(Matrix)
library(stats)
library(cjoint)
library(estimatr)
library(conjoint)
library(dotwhisker)
library(sjPlot)
library(knitr)
library(cregg)

#IMPORTANT NOTE: may need to install cregg from archive
#url <- "https://cran.r-project.org/src/contrib/Archive/cregg/cregg_0.4.0.tar.gz"
#pkgFile <- "cregg_0.4.0.tar.gz"
#download.file(url = url, destfile = pkgFile)
# expand zip file "cregg"
# look at the DESCRIPTION file in the expanded package directory
# Install dependencies list in the DESCRIPTION file
#install.packages(pkgs=pkgFile, type="source", repos=NULL)
#unlink(pkgFile)
#if desired: delete zip file and subsequential file

##
##
##



# load primary data

LiberiaRound1And2 <- read_excel("/Users/cmailhot/Desktop/Replication of Mailhot and Karim. International Statebuilding and Civilian Preferences. IO/Data.xlsx")

# initial data clean-up 

#drop respondents who didn't get to conjoint:

LiberiaRound1And2 <- subset(LiberiaRound1And2, ActorAssignedNumber !="NA")

#create numerical representations of the levels of each attribute (actor, action, institution)
LiberiaRound1And2$Actor[LiberiaRound1And2$ActorAssignedNumber == "Liberian Gov't"] <- 1
LiberiaRound1And2$Actor[LiberiaRound1And2$ActorAssignedNumber == "United Nations"] <- 2
LiberiaRound1And2$Actor[LiberiaRound1And2$ActorAssignedNumber == "United States"] <- 3
LiberiaRound1And2$Actor[LiberiaRound1And2$ActorAssignedNumber == "China"] <- 4

#convert to factor variable
LiberiaRound1And2$Actor <- factor(LiberiaRound1And2$Actor, levels = c(1,2,3,4), labels = c("Liberian Gov't", "United Nations", "United States", "China"))
is.factor(LiberiaRound1And2$Actor)

LiberiaRound1And2$Action[LiberiaRound1And2$ActionAssignedNumber == "Rebuilt"] <- 1
LiberiaRound1And2$Action[LiberiaRound1And2$ActionAssignedNumber == "Monitored Reconstruction"] <- 2

#convert to factor variable
LiberiaRound1And2$Action <- factor(LiberiaRound1And2$Action, levels = c(1,2), labels = c("Rebuilt", "Monitored Reconstruction"))
is.factor(LiberiaRound1And2$Action)


LiberiaRound1And2$Institution[LiberiaRound1And2$InstitutionAssignedNumber == "Postal Service"] <- 1
LiberiaRound1And2$Institution[LiberiaRound1And2$InstitutionAssignedNumber == "Land Authority"] <- 1
LiberiaRound1And2$Institution[LiberiaRound1And2$InstitutionAssignedNumber == "AFL Special Forces Team"] <- 2

#convert to factor variable
LiberiaRound1And2$Institution <- factor(LiberiaRound1And2$Institution, levels = c(1,2), labels = c("Postal Service/Land Authority", "AFL Special Forces Team"))
is.factor(LiberiaRound1And2$Institution)

#
#
#
# Analyses in Main Manuscript
#
#
#

model <- ComboPreference ~ Actor + Action + Institution 

# general results
generalmarginalmeanslabels <- c("Postal Service/Land Authority", "AFL Special Forces Team", "INSTITUTION", "Rebuilt", "Monitored Reconstruction", "ACTION", "Liberian Gov't", "United Nations", "United States", "China", "ACTOR")
generalmarginalmeans <- cj(LiberiaRound1And2, model, id = ~result_id, estimate = "mm")

# Marginal means: Figure 1 (main manuscript) and Figure A4 (appendix)
plot(generalmarginalmeans, vline = 0.5) + labs(color = "Attributes", x = "Marginal Means", y  = "Statebuilding Characteristics") + scale_colour_grey(start = 0, end = .75, na.translate = FALSE) +  theme_bw()  +  theme(legend.position = "bottom", legend.title = element_blank()) + scale_x_continuous(minor_breaks = seq(0.25, 0.75, .25)) + scale_y_discrete(labels = generalmarginalmeanslabels, limits =c("Postal Service/Land Authority", "AFL Special Forces Team", "INSTITUTION", "Rebuilt", "Monitored Reconstruction", "ACTION", "Liberian Gov't", "United Nations", "United States", "China", "ACTOR")) 


# Marginal Means: Table A1 (appendix)
head(generalmarginalmeans[c("feature", "level", "estimate", "std.error")], 8L)
generalmarginalmeans <- subset(generalmarginalmeans, select = c("feature", "level", "estimate", "std.error", "z", "p"))
kable(generalmarginalmeans, "latex", booktabs = TRUE)


#AMCEs: Table A5 (appendix)
generalamces <- cj(LiberiaRound1And2, model, id = ~result_id, estimate="amce")
generalamces <- subset(generalamces, select = c("feature", "level", "estimate", "std.error", "z", "p"))
kable(generalamces, "latex", booktabs = TRUE)


#
#Results by sex
#

LiberiaRound1And2$Sex <- factor(LiberiaRound1And2$back_sex_Male, levels = c(0, 1), labels = c("Female", "Male"))
generalmarginalmeans_bysex <- cj(LiberiaRound1And2, model, id = ~result_id, estimate = "mm", by = ~Sex)

#Marginal Means: Figure 2, part 1 (main manuscript)
#note: pink for female; blue for male
plot(generalmarginalmeans_bysex, vline = 0.5, group = "Sex", xlim=c(0.4,0.6)) + ggtitle("Sex") + labs(color = "Attributes", x = "Marginal Means", y  = "Attributes and Levels") + scale_colour_grey(start = 0, end = .75, na.translate = FALSE) +  theme_bw()  +  theme(legend.position = "bottom", plot.title = element_text(hjust = 0.49)) + scale_y_discrete(labels = generalmarginalmeanslabels, limits =c("Postal Service/Land Authority", "AFL Special Forces Team", "INSTITUTION", "Rebuilt", "Monitored Reconstruction", "ACTION", "Liberian Gov't", "United Nations", "United States", "China", "ACTOR")) 

#Marginal Means: Table A2 (appendix)
#Note: First set of results in table for female respondents, second set of results in table for male respondents
generalmarginalmeans_bysex <- subset(generalmarginalmeans_bysex, select = c("feature", "level", "estimate", "std.error", "z", "p"))
kable(generalmarginalmeans_bysex, "latex", booktabs = TRUE)

# AMCEs: Table A6 (appendix)
#Note: First set of results in table for female respondents, second set of results in table for male respondents
generalamces_bysex <- cj(LiberiaRound1And2, model, id = ~result_id, estimate = "amce", by = ~Sex)
generalamces_bysex <- subset(generalamces_bysex, select = c("feature", "level", "estimate", "std.error", "z", "p"))
kable(generalamces_bysex, "latex", booktabs = TRUE)


#
#Results by Income
#
LiberiaRound1And2Class <- subset(LiberiaRound1And2, back_income <= 30)
LiberiaRound1And2Class$Class <- ifelse(LiberiaRound1And2Class$back_income <= 5, "Low Income", "High Income")
LiberiaRound1And2Class$Class <- factor(LiberiaRound1And2Class$Class)
LiberiaRound1And2Class$Class <- factor(LiberiaRound1And2Class$Class)

generalmarginalmeans_byclass <- cj(LiberiaRound1And2Class, model, id = ~result_id, estimate = "mm", by = ~Class)

#Marginal Means: Figure 2, part 2 (main manuscript)
#note: blue: low income; pink: high income
plot(generalmarginalmeans_byclass, vline = 0.5, group = "Class", xlim=c(0.4,0.6)) + ggtitle("Class") + labs(color = "Attributes", x = "Marginal Means", y  = "Attributes and Levels") +  scale_colour_grey(start = 0, end = .75, na.translate = FALSE) +  theme_bw()  +  theme(legend.position = "bottom", plot.title = element_text(hjust = 0.49)) + scale_y_discrete(labels = generalmarginalmeanslabels, limits =c("Postal Service/Land Authority", "AFL Special Forces Team", "INSTITUTION", "Rebuilt", "Monitored Reconstruction", "ACTION", "Liberian Gov't", "United Nations", "United States", "China", "ACTOR")) 

#Marginal Means: Table A3 (appendix)
#Note: First set of results in table for high income respondents, second set of results in table for low income respondents
generalmarginalmeans_byclass <- subset(generalmarginalmeans_byclass, select = c("feature", "level", "estimate", "std.error", "z", "p"))
kable(generalmarginalmeans_byclass, "latex", booktabs = TRUE)

# AMCEs: Table A7 (appendix)
#Note: First set of results in table for high income respondents, second set of results in table for low income respondents
generalamces_byclass <- cj(LiberiaRound1And2Class, model, id = ~result_id, estimate = "amce", by = ~Class)
generalamces_byclass <- subset(generalamces_byclass, select = c("feature", "level", "estimate", "std.error", "z", "p"))
kable(generalamces_byclass, "latex", booktabs = TRUE)

#
#Results by Wartime Conflict Experiences
#

Round1 <- subset(LiberiaRound1And2, Round1_2020 == 1)
Round1$Actor <- 0
Round1$Action <- 0
Round1$Institution <- 0
#create numerical representations of the levels of each attribute (actor, action, institution)
Round1$Actor[Round1$ActorAssignedNumber == "Liberian Gov't"] <- 1
Round1$Actor[Round1$ActorAssignedNumber == "United Nations"] <- 2
Round1$Actor[Round1$ActorAssignedNumber == "United States"] <- 3
Round1$Actor[Round1$ActorAssignedNumber == "China"] <- 4

#convert to factor variable
Round1$Actor <- factor(Round1$Actor, levels = c(1,2,3,4), labels = c("Liberian Gov't", "United Nations", "United States", "China"))
is.factor(Round1$Actor)

Round1$Action[Round1$ActionAssignedNumber == "Rebuilt"] <- 1
Round1$Action[Round1$ActionAssignedNumber == "Monitored Reconstruction"] <- 2

#convert to factor variable
Round1$Action <- factor(Round1$Action, levels = c(1,2), labels = c("Rebuilt", "Monitored Reconstruction"))
is.factor(Round1$Action)


Round1$Institution[Round1$InstitutionAssignedNumber == "Postal Service"] <- 1
Round1$Institution[Round1$InstitutionAssignedNumber == "AFL Special Forces Team"] <- 2

#convert to factor variable
Round1$Institution <- factor(Round1$Institution, levels = c(1,2), labels = c("Postal Service", "AFL Special Forces Team"))
is.factor(Round1$Institution)

#create variable that defines levels of wartime community violence experience from 3 other variables 

Round1['ConflictExposureLevels'] <- NA
Round1$ConflictExposureLevels[Round1$back_neighbhdbad_NotAffectedAtAl ==1] = "Not Affected"
Round1$ConflictExposureLevels[Round1$back_neighbhdbad_AffectedToMildD ==1] = "Mild Degree"
Round1$ConflictExposureLevels[Round1$back_neighbhdbad_AffectedToModer ==1] = "Moderate Degree"
Round1$ConflictExposureLevels[Round1$back_neighbhdbad_AffectedToSever ==1] = "Severe Degree"

LiberiaRound1ConflictExposure <- subset(Round1, ConflictExposureLevels != "NA")
LiberiaRound1ConflictExposure$ConflictExposureLevels <- factor(LiberiaRound1ConflictExposure$ConflictExposureLevels, levels = c("Not Affected", "Mild Degree", "Moderate Degree", "Severe Degree"))

generalmarginalmeanslabelsround1 <- c("Postal Office", "AFL Special Forces Team", "INSTITUTION", "Rebuilt", "Monitored Reconstruction", "ACTION", "Liberian Gov't", "United Nations", "United States", "China", "ACTOR")
generalmarginalmeans_bywartimeexposure <- cj(LiberiaRound1ConflictExposure, model, id = ~result_id, estimate = "mm", by = ~ConflictExposureLevels)

#Marginal Means: Figure 2, part 3 (main manuscript)
#note: results for "not affected" in orange; "mild degree" in green, "moderate degree" in blue, "severe degree" in purple"
plot(generalmarginalmeans_bywartimeexposure, vline = 0.5, group = "ConflictExposureLevels", xlim=c(0.2, 0.8)) + ggtitle("Wartime Community Impacts") + labs(color = "Attributes", x = "Marginal Means", y  = "Attributes and Levels") + scale_colour_grey(start = 0, end = .75, na.translate = FALSE) +  theme_bw()  +  theme(legend.title = element_blank(), legend.text = element_text(size = 7), legend.position = "bottom", plot.title = element_text(hjust = 0.49)) + scale_y_discrete(labels = generalmarginalmeanslabelsround1, limits =c("Postal Service", "AFL Special Forces Team", "INSTITUTION", "Rebuilt", "Monitored Reconstruction", "ACTION", "Liberian Gov't", "United Nations", "United States", "China", "ACTOR")) 


#Marginal Means: Table A4 (appendix)
#note: First set of results: not affected; second set of results: mild degree; third set of results: moderate degree; fourth set of results: severe degree
generalmarginalmeans_bywartimeexposure <- subset(generalmarginalmeans_bywartimeexposure, select = c("feature", "level", "estimate", "std.error", "z", "p"))
kable(generalmarginalmeans_bywartimeexposure, "latex", booktabs = TRUE)


#AMCEs: Table A8 (appendix)
generalamces_bywartimeexperiences <- cj(LiberiaRound1ConflictExposure, model, id = ~result_id, estimate = "amce", by = ~ConflictExposureLevels)
generalamces_bywartimeexperiences <- subset(generalamces_bywartimeexperiences, select = c("feature", "level", "estimate", "std.error", "z", "p"))
kable(generalamces_bywartimeexperiences, "latex", booktabs = TRUE)



#
#
# Analyses only in Appendix
#
#
#

# Location 

LiberiaRound1And2$Location <- factor(LiberiaRound1And2$town_name, levels = c(1, 2, 3), labels = c("Monrovia", "Gbarnga", "Bong Mines"))

#AMCEs: Table A9 (appendix)
#Note: First set of results for Monrovia, second set of results for Gbarnga, third set of results for Bong Mines
generalamces_bylocation <- cj(LiberiaRound1And2, model, id = ~result_id, estimate = "amce", by = ~Location)
generalamces_bylocation <- subset(generalamces_bylocation, select = c("feature", "level", "estimate", "std.error", "z", "p"))
kable(generalamces_bylocation, "latex", booktabs = TRUE)


# education 


LiberiaRound1And2['Education'] <- NA
LiberiaRound1And2$Education[LiberiaRound1And2$back_education <5] = "Primary Education"
LiberiaRound1And2$Education[LiberiaRound1And2$back_education == 5] = "High School Education"
LiberiaRound1And2$Education[LiberiaRound1And2$back_education > 5] = "University/Higher Education"
LiberiaRound1And2$Education <- factor(LiberiaRound1And2$Education)
#Marginal Means: Table A10 (appendix)
#Note: First set of results for high school, second set of results for primary school, third set of results for university education 
generalmarginalmeans_byeducation <- cj(LiberiaRound1And2, model, id = ~result_id, estimate = "mm", by = ~Education)
generalmarginalmeans_byeducation <- subset(generalmarginalmeans_byeducation, select = c("feature", "level", "estimate", "std.error", "z", "p"))
kable(generalmarginalmeans_byeducation, "latex", booktabs = TRUE)



# Location 

LiberiaRound1And2$Location <- factor(LiberiaRound1And2$town_name, levels = c(1, 2, 3), labels = c("Monrovia", "Gbarnga", "Bong Mines"))

#Marginal Means: Table A11 (appendix)
#Note: First set of results for Monrovia, second set of results for Gbarnga, third set of results for Bong Mines
generalmarginalmeans_bylocation <- cj(LiberiaRound1And2, model, id = ~result_id, estimate = "mm", by = ~Location)
generalmarginalmeans_bylocation <- subset(generalmarginalmeans_bylocation, select = c("feature", "level", "estimate", "std.error", "z", "p"))
kable(generalmarginalmeans_bylocation, "latex", booktabs = TRUE)


#ethnicity 
LiberiaRound1And2['Ethnicity'] <- NA
LiberiaRound1And2$Ethnicity[LiberiaRound1And2$back_ethnicity == 1] = "Bassa"
LiberiaRound1And2$Ethnicity[LiberiaRound1And2$back_ethnicity == 10] = "Kpelle"
LiberiaRound1And2$Ethnicity[LiberiaRound1And2$back_ethnicity == 14] = "Mandingo"
LiberiaRound1And2$Ethnicity[LiberiaRound1And2$back_ethnicity == 2 | LiberiaRound1And2$back_ethnicity == 3 | LiberiaRound1And2$back_ethnicity == 4 | LiberiaRound1And2$back_ethnicity == 5 | LiberiaRound1And2$back_ethnicity == 6 | LiberiaRound1And2$back_ethnicity == 7 | LiberiaRound1And2$back_ethnicity == 8 | LiberiaRound1And2$back_ethnicity == 9 | LiberiaRound1And2$back_ethnicity == 11 | LiberiaRound1And2$back_ethnicity == 12 | LiberiaRound1And2$back_ethnicity == 13 | LiberiaRound1And2$back_ethnicity == 15 | LiberiaRound1And2$back_ethnicity == 16 | LiberiaRound1And2$back_ethnicity == 17 ] = "Other"

LiberiaRound1And2$Ethnicity <- factor(LiberiaRound1And2$Ethnicity)

#Marginal Means: Table A12 (appendix)
#Note: First set of results for Bassa, second set of results for Kpelle, third set of results for mandingo, fourth set of results for "Other" 
generalmarginalmeans_byethnicity <- cj(LiberiaRound1And2, model, id = ~result_id, estimate = "mm", by = ~Ethnicity)
generalmarginalmeans_byethnicity <- subset(generalmarginalmeans_byethnicity, select = c("feature", "level", "estimate", "std.error", "z", "p"))
kable(generalmarginalmeans_byethnicity, "latex", booktabs = TRUE)


# location type (Urban/Rural)


# urban/rural divide

LiberiaRound1And2$Capital1Rural2 <- factor(LiberiaRound1And2$Capital1Rural2, levels = c(1,2), labels = c("Rural", "Urban (Capital)"))
is.factor(LiberiaRound1And2$Capital1Rural2)

#Marginal Means: Table A13 (appendix)
#Note: First set of resuls for rural, second set of results for urban
generalmarginalmeans_byurbanrural <- cj(LiberiaRound1And2, model, id = ~result_id, estimate = "mm", by = ~Capital1Rural2)
generalmarginalmeans_byurbanrural <- subset(generalmarginalmeans_byurbanrural, select = c("feature", "level", "estimate", "std.error", "z", "p"))
kable(generalmarginalmeans_byurbanrural, "latex", booktabs = TRUE)


# Survey Round

LiberiaRound1And2$Round <- 1
LiberiaRound1And2$Round[LiberiaRound1And2$Round2_2022 == 1] <- 2
LiberiaRound1And2$Round <- factor(LiberiaRound1And2$Round, levels = c(1,2), labels = c("Round 1", "Round 2"))

#Marginal Means: Table A14 (appendix)
#Note: First set of resuls for Round 1, second set of results for Round 2
generalmarginalmeans_byround <- cj(LiberiaRound1And2, model, id = ~result_id, estimate = "mm", by = ~Round)
generalmarginalmeans_byround <- subset(generalmarginalmeans_byround, select = c("feature", "level", "estimate", "std.error", "z", "p"))
kable(generalmarginalmeans_byround, "latex", booktabs = TRUE)



# Political Party 

LiberiaRound1And2['PoliticalParty'] <- NA
LiberiaRound1And2$PoliticalParty[LiberiaRound1And2$pol_party == 1 | LiberiaRound1And2$pol_party == 2 | LiberiaRound1And2$pol_party == 3 | LiberiaRound1And2$pol_party == 4 | LiberiaRound1And2$pol_party == 5 | LiberiaRound1And2$pol_party == 6 | LiberiaRound1And2$pol_party == 7 | LiberiaRound1And2$pol_party == 8 | LiberiaRound1And2$pol_party == 9 | LiberiaRound1And2$pol_party == 10 | LiberiaRound1And2$pol_party == 14 | LiberiaRound1And2$pol_party == 98 | LiberiaRound1And2$pol_party == 99] = "Other"
LiberiaRound1And2$PoliticalParty[LiberiaRound1And2$pol_party == 11] = "Unity Party"
LiberiaRound1And2$PoliticalParty[LiberiaRound1And2$pol_party == 12] = "CDC"
LiberiaRound1And2$PoliticalParty <- factor(LiberiaRound1And2$PoliticalParty)


#Marginal Means: Table A15 (appendix)
#Note: First set of resuls for CDC , second set of results for "Other," third set of results for Unity Party
generalmarginalmeans_bypoliticalparty <- cj(LiberiaRound1And2, model, id = ~result_id, estimate = "mm", by = ~PoliticalParty)
generalmarginalmeans_bypoliticalparty <- subset(generalmarginalmeans_bypoliticalparty, select = c("feature", "level", "estimate", "std.error", "z", "p"))
kable(generalmarginalmeans_bypoliticalparty, "latex", booktabs = TRUE)



# Power Analyses 

df2 = cjpowr_amce(amce = 0.05, power = 0.8, levels = 4)

df2$n/(2*4)

cjpowr_amce(amce = 0.05, n = 9780, levels = 4)
# --> Power of 0.94, greater than 0.80 (conventional)


cjpowr_amce_vec <- Vectorize(cjpowr_amce)
d <- expand.grid(amce = c(0.05), n = seq(from = 100, to = 50000, length.out = 1000))

df <- t(cjpowr_amce_vec(amce = d$amce, n = d$n, sims = 100000, levels = 4, alpha = 0.05))
df <- data.frame(df)
df[] <- lapply(df, unlist)

library(plotly)


#Figure A1 and A3 (appendix)
plot_ly(df, x = ~n, y = ~power, type = 'scatter', mode = 'lines', linetype = ~amce) %>%
  layout(
    xaxis = list(title = "Minimum Effective Sample Size",
                 zeroline = F,
                 hoverformat = '.0f'),
    yaxis = list(title = "Power",
                 range = c(0,1),
                 zeroline = F,
                 hoverformat = '.2f'),
    legend=list(title=list(text='<b> AMCE </b>')),
    hovermode = "x unified"
  )



# balance tests

LiberiaRound1And2$China = 0 
LiberiaRound1And2$China[LiberiaRound1And2$Actor == "China"] <- 1
LiberiaRound1And2$LiberianGovt = 0 
LiberiaRound1And2$LiberianGovt[LiberiaRound1And2$Actor == "Liberian Gov't"] <- 1
LiberiaRound1And2$UnitedNations = 0 
LiberiaRound1And2$UnitedNations[LiberiaRound1And2$Actor == "United Nations"] <- 1
LiberiaRound1And2$UnitedStates = 0 
LiberiaRound1And2$UnitedStates[LiberiaRound1And2$Actor == "United States"] <- 1
LiberiaRound1And2$Rebuilt = 0
LiberiaRound1And2$Rebuilt[LiberiaRound1And2$Action == "Rebuilt"] <- 1
LiberiaRound1And2$MonitoredReconstruction = 0
LiberiaRound1And2$MonitoredReconstruction[LiberiaRound1And2$Action == "Monitored Reconstruction"] <- 1
LiberiaRound1And2$AFLSpecialForces = 0
LiberiaRound1And2$AFLSpecialForces[LiberiaRound1And2$Institution == "AFL Special Forces Team"] <- 1
LiberiaRound1And2$PostalServiceLandAuthority = 0
LiberiaRound1And2$PostalServiceLandAuthority[LiberiaRound1And2$Institution == "Postal Service/Land Authority"] <- 1


LiberiaRound1And2Class$China = 0 
LiberiaRound1And2Class$China[LiberiaRound1And2Class$Actor == "China"] <- 1
LiberiaRound1And2Class$LiberianGovt = 0 
LiberiaRound1And2Class$LiberianGovt[LiberiaRound1And2Class$Actor == "Liberian Gov't"] <- 1
LiberiaRound1And2Class$UnitedNations = 0 
LiberiaRound1And2Class$UnitedNations[LiberiaRound1And2Class$Actor == "United Nations"] <- 1
LiberiaRound1And2Class$UnitedStates = 0 
LiberiaRound1And2Class$UnitedStates[LiberiaRound1And2Class$Actor == "United States"] <- 1
LiberiaRound1And2Class$Rebuilt = 0
LiberiaRound1And2Class$Rebuilt[LiberiaRound1And2Class$Action == "Rebuilt"] <- 1
LiberiaRound1And2Class$MonitoredReconstruction = 0
LiberiaRound1And2Class$MonitoredReconstruction[LiberiaRound1And2Class$Action == "Monitored Reconstruction"] <- 1
LiberiaRound1And2Class$AFLSpecialForces = 0
LiberiaRound1And2Class$AFLSpecialForces[LiberiaRound1And2Class$Institution == "AFL Special Forces Team"] <- 1
LiberiaRound1And2Class$PostalServiceLandAuthority = 0
LiberiaRound1And2Class$PostalServiceLandAuthority[LiberiaRound1And2Class$Institution == "Postal Service/Land Authority"] <- 1

LiberiaRound1ConflictExposure$China = 0 
LiberiaRound1ConflictExposure$China[LiberiaRound1ConflictExposure$Actor == "China"] <- 1
LiberiaRound1ConflictExposure$LiberianGovt = 0 
LiberiaRound1ConflictExposure$LiberianGovt[LiberiaRound1ConflictExposure$Actor == "Liberian Gov't"] <- 1
LiberiaRound1ConflictExposure$UnitedNations = 0 
LiberiaRound1ConflictExposure$UnitedNations[LiberiaRound1ConflictExposure$Actor == "United Nations"] <- 1
LiberiaRound1ConflictExposure$UnitedStates = 0 
LiberiaRound1ConflictExposure$UnitedStates[LiberiaRound1ConflictExposure$Actor == "United States"] <- 1
LiberiaRound1ConflictExposure$Rebuilt = 0
LiberiaRound1ConflictExposure$Rebuilt[LiberiaRound1ConflictExposure$Action == "Rebuilt"] <- 1
LiberiaRound1ConflictExposure$MonitoredReconstruction = 0
LiberiaRound1ConflictExposure$MonitoredReconstruction[LiberiaRound1ConflictExposure$Action == "Monitored Reconstruction"] <- 1
LiberiaRound1ConflictExposure$AFLSpecialForces = 0
LiberiaRound1ConflictExposure$AFLSpecialForces[LiberiaRound1ConflictExposure$Institution == "AFL Special Forces Team"] <- 1
LiberiaRound1ConflictExposure$PostalService = 0
LiberiaRound1ConflictExposure$PostalService[LiberiaRound1ConflictExposure$Institution == "Postal Service"] <- 1


### Table A16 (appendix)

LiberiaRound1And2 %>% group_by(Sex) %>% summarize(mean = mean(China), sum = sum(China), sd = sd(China))
LiberiaRound1And2 %>% group_by(Sex) %>% summarize(mean = mean(LiberianGovt), sum = sum(LiberianGovt), sd = sd(LiberianGovt))
LiberiaRound1And2 %>% group_by(Sex) %>% summarize(mean = mean(UnitedNations), sum = sum(UnitedNations), sd = sd(UnitedNations))
LiberiaRound1And2 %>% group_by(Sex) %>% summarize(mean = mean(UnitedStates), sum = sum(UnitedStates), sd = sd(UnitedStates))
LiberiaRound1And2 %>% group_by(Sex) %>% summarize(mean = mean(Rebuilt), sum = sum(Rebuilt), sd = sd(Rebuilt))
LiberiaRound1And2 %>% group_by(Sex) %>% summarize(mean = mean(MonitoredReconstruction), sum = sum(MonitoredReconstruction), sd = sd(MonitoredReconstruction))
LiberiaRound1And2 %>% group_by(Sex) %>% summarize(mean = mean(AFLSpecialForces), sum = sum(AFLSpecialForces), sd = sd(AFLSpecialForces))
LiberiaRound1And2 %>% group_by(Sex) %>% summarize(mean = mean(PostalServiceLandAuthority), sum = sum(PostalServiceLandAuthority), sd = sd(PostalServiceLandAuthority))

LiberiaRound1And2Class %>% group_by(Class) %>% summarize(mean = mean(China), sum = sum(China), sd = sd(China))
LiberiaRound1And2Class %>% group_by(Class) %>% summarize(mean = mean(LiberianGovt), sum = sum(LiberianGovt), sd = sd(LiberianGovt))
LiberiaRound1And2Class %>% group_by(Class) %>% summarize(mean = mean(UnitedNations), sum = sum(UnitedNations), sd = sd(UnitedNations))
LiberiaRound1And2Class %>% group_by(Class) %>% summarize(mean = mean(UnitedStates), sum = sum(UnitedStates), sd = sd(UnitedStates))
LiberiaRound1And2Class %>% group_by(Class) %>% summarize(mean = mean(Rebuilt), sum = sum(Rebuilt), sd = sd(Rebuilt))
LiberiaRound1And2Class %>% group_by(Class) %>% summarize(mean = mean(MonitoredReconstruction), sum = sum(MonitoredReconstruction), sd = sd(MonitoredReconstruction))
LiberiaRound1And2Class %>% group_by(Class) %>% summarize(mean = mean(AFLSpecialForces), sum = sum(AFLSpecialForces), sd = sd(AFLSpecialForces))
LiberiaRound1And2Class %>% group_by(Class) %>% summarize(mean = mean(PostalServiceLandAuthority), sum = sum(PostalServiceLandAuthority), sd = sd(PostalServiceLandAuthority))

LiberiaRound1ConflictExposure %>% group_by(ConflictExposureLevels) %>% summarize(mean = mean(China), sum = sum(China), sd = sd(China))
LiberiaRound1ConflictExposure %>% group_by(ConflictExposureLevels) %>% summarize(mean = mean(LiberianGovt), sum = sum(LiberianGovt), sd = sd(LiberianGovt))
LiberiaRound1ConflictExposure %>% group_by(ConflictExposureLevels) %>% summarize(mean = mean(UnitedNations), sum = sum(UnitedNations), sd = sd(UnitedNations))
LiberiaRound1ConflictExposure %>% group_by(ConflictExposureLevels) %>% summarize(mean = mean(UnitedStates), sum = sum(UnitedStates), sd = sd(UnitedStates))
LiberiaRound1ConflictExposure %>% group_by(ConflictExposureLevels) %>% summarize(mean = mean(Rebuilt), sum = sum(Rebuilt), sd = sd(Rebuilt))
LiberiaRound1ConflictExposure %>% group_by(ConflictExposureLevels) %>% summarize(mean = mean(MonitoredReconstruction), sum = sum(MonitoredReconstruction), sd = sd(MonitoredReconstruction))
LiberiaRound1ConflictExposure %>% group_by(ConflictExposureLevels) %>% summarize(mean = mean(AFLSpecialForces), sum = sum(AFLSpecialForces), sd = sd(AFLSpecialForces))
LiberiaRound1ConflictExposure %>% group_by(ConflictExposureLevels) %>% summarize(mean = mean(PostalService), sum = sum(PostalService), sd = sd(PostalService))

LiberiaRound1And2 %>% group_by(Education) %>% summarize(mean = mean(China), sum = sum(China), sd = sd(China))
LiberiaRound1And2 %>% group_by(Education) %>% summarize(mean = mean(LiberianGovt), sum = sum(LiberianGovt), sd = sd(LiberianGovt))
LiberiaRound1And2 %>% group_by(Education) %>% summarize(mean = mean(UnitedNations), sum = sum(UnitedNations), sd = sd(UnitedNations))
LiberiaRound1And2 %>% group_by(Education) %>% summarize(mean = mean(UnitedStates), sum = sum(UnitedStates), sd = sd(UnitedStates))
LiberiaRound1And2 %>% group_by(Education) %>% summarize(mean = mean(Rebuilt), sum = sum(Rebuilt), sd = sd(Rebuilt))
LiberiaRound1And2 %>% group_by(Education) %>% summarize(mean = mean(MonitoredReconstruction), sum = sum(MonitoredReconstruction), sd = sd(MonitoredReconstruction))
LiberiaRound1And2 %>% group_by(Education) %>% summarize(mean = mean(AFLSpecialForces), sum = sum(AFLSpecialForces), sd = sd(AFLSpecialForces))
LiberiaRound1And2 %>% group_by(Education) %>% summarize(mean = mean(PostalServiceLandAuthority), sum = sum(PostalServiceLandAuthority), sd = sd(PostalServiceLandAuthority))

LiberiaRound1And2 %>% group_by(Ethnicity) %>% summarize(mean = mean(China), sum = sum(China), sd = sd(China))
LiberiaRound1And2 %>% group_by(Ethnicity) %>% summarize(mean = mean(LiberianGovt), sum = sum(LiberianGovt), sd = sd(LiberianGovt))
LiberiaRound1And2 %>% group_by(Ethnicity) %>% summarize(mean = mean(UnitedNations), sum = sum(UnitedNations), sd = sd(UnitedNations))
LiberiaRound1And2 %>% group_by(Ethnicity) %>% summarize(mean = mean(UnitedStates), sum = sum(UnitedStates), sd = sd(UnitedStates))
LiberiaRound1And2 %>% group_by(Ethnicity) %>% summarize(mean = mean(Rebuilt), sum = sum(Rebuilt), sd = sd(Rebuilt))
LiberiaRound1And2 %>% group_by(Ethnicity) %>% summarize(mean = mean(MonitoredReconstruction), sum = sum(MonitoredReconstruction), sd = sd(MonitoredReconstruction))
LiberiaRound1And2 %>% group_by(Ethnicity) %>% summarize(mean = mean(AFLSpecialForces), sum = sum(AFLSpecialForces), sd = sd(AFLSpecialForces))
LiberiaRound1And2 %>% group_by(Ethnicity) %>% summarize(mean = mean(PostalServiceLandAuthority), sum = sum(PostalServiceLandAuthority), sd = sd(PostalServiceLandAuthority))


